/* Declaration of set of functions to solve numerictly */
#include "jkFunc.h"
using namespace std;

/* FROM JK AND JJ PAPER - SGA METOD */
// In Par parameter:
//  0. h (external magnetic field),
//  1. n (band filling),
//  2. t1 (hopping integral for the nearest neighbors),
//  3. t2 (hopping integral for the second nearest neighbors),
//  4. beta (inverse temperature),
//  5. U (Hubbard parameter [repulsive potential for electrons in the same
//     spatial point]).
//  6. L (size of lattice)
//
// In x list:
//  0. mu (chemical potential)
//  1. m (magnetization)
//  2. d (density of double occupied sities)
//  3. lambda_n (for SGA lagrange multiplier connected with n)
//  4. lambda_m (for SGA lagrange multiplier connected with m)


int jk2dGA (const gsl_vector *_x, void *_Par, gsl_vector * fun)
{
	Parameters *Par = reinterpret_cast<Parameters*>(_Par);

	const unsigned int N = 3; // we have set of 3 functions to solve
	double x[N]; // vector of arguments
	for(unsigned int i=0; i<N; i++) { x[i] = gsl_vector_get(_x, i); }
	JK *jk = new JKGA(Par);
	double *y = jk->sumation2D(x);

	// I put finally result from y to gsl vector "fun"
   for(unsigned int i=0; i<N; i++) { gsl_vector_set (fun, i, y[i]); }

   delete jk;
   return GSL_SUCCESS;
}

int jk2dSGA (const gsl_vector *_x, void *_Par, gsl_vector * fun)
{
	Parameters *Par = reinterpret_cast<Parameters*>(_Par);

	const unsigned int N = 5; // we have set of 3 functions to solve
	double x[N]; // vector of arguments
	for(unsigned int i=0; i<N; i++) { x[i] = gsl_vector_get(_x, i); }
	JK *jk = new JKSGA(Par);
	double *y = jk->sumation2D(x);

//	printf("x = ");
//	for(unsigned int i=0; i<N; i++) printf("%e ", x[i]);
//	printf("y = ");
//	for(unsigned int i=0; i<N; i++) printf("%e ", y[i]);
//	printf("\n");

	// I put finally result from y to gsl vector "fun"
   for(unsigned int i=0; i<N; i++) { gsl_vector_set (fun, i, y[i]); }

   delete jk;
   return GSL_SUCCESS;
}

int jk3dGA (const gsl_vector *_x, void *_Par, gsl_vector * fun)
{
	Parameters *Par = reinterpret_cast<Parameters*>(_Par);

	const unsigned int N = 3; // we have set of 3 functions to solve
	double x[N]; // vector of arguments
	for(unsigned int i=0; i<N; i++) { x[i] = gsl_vector_get(_x, i); }
	JK *jk = new JKGA(Par);
	double *y = jk->sumation3D(x);

	// I put finally result from y to gsl vector "fun"
   for(unsigned int i=0; i<N; i++) { gsl_vector_set (fun, i, y[i]); }

   delete jk;
   return GSL_SUCCESS;
}

int jk3dSGA (const gsl_vector *_x, void *_Par, gsl_vector * fun)
{
	Parameters *Par = reinterpret_cast<Parameters*>(_Par);

	const unsigned int N = 5; // we have set of 3 functions to solve
	double x[N]; // vector of arguments
	for(unsigned int i=0; i<N; i++) { x[i] = gsl_vector_get(_x, i); }
	JK *jk = new JKSGA(Par);
	double *y = jk->sumation3D(x);

	// I put finally result from y to gsl vector "fun"
   for(unsigned int i=0; i<N; i++) { gsl_vector_set (fun, i, y[i]); }

   delete jk;
   return GSL_SUCCESS;
}

/***********************************************************************/

